Mean-field dynamics of a non-hermitian Bose-Hubbard dimer 
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We investigate an Af-particle Bose-Hubbard dimer with an additional effective decay term in one of the sites. 
A mean-field approximation for this non-hermitian many-particle system is derived, based on a coherent state 
approximation. The resulting nonlinear, non-hermitian two-level dynamics, in particular the fixed point struc- 
tures showing characteristic modifications of the self trapping transition, are analyzed. The mean-field dynamics 
is found to be in reasonable agreement with the full many-particle evolution. 
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In the theoretical investigation of Bose-Einstein conden- 
sates (BEC) the mean-field approximation leading to the de- 
scription via a Gross-Pitaevskii nonlinear Schrodinger equa- 
tion (GPE) is almost indispensable. It is usually achieved 
by replacing the bosonic field operators in the multi-particle 
system with c-numbers (the effective single-particle conden- 
sate wave functions), and describes the system quite well for 
large particle numbers and low temperatures. This approach 
is closely related to a classicalization [1] and allows for the 
application of semiclassical methods [2] . 

Recently considerable attention has been paid to effective 
non-hermitian mean-field theories describing the scattering 
and transport behavior of BECs Ok as well as the implications 
of decay (boundary dissipation) |4J,|5|,|6|]. The latter is closely 
related to atom laser, for which it is possible to go beyond the 
mean-field approximation and calculate the eigenmodes us- 
ing Fano diagonalization [7]. For linear quantum systems, an 
effective non-hermitian Hamiltonian formalism proved useful 
and instructive for the description of open quantum systems 
in various fields of physics. Non-hermitian Hamiltonians typ- 
ically yield complex eigenvalues whose imaginary parts de- 
scribe the rates with which an eigenstate decays to the external 
world. Other kinds of non-hermitian (PT- symmetric) quantum 
theories have also been suggested as a generalization of quan- 
tum mechanics on a fundamental level |8|]. 

However, the non-hermitian GPE has been formulated in an 
ad hoc manner as a generalization of the mean-field Hamil- 
tonian and a derivation starting from a non-hermitian many 
particle system is required. This is as well interesting in a 
wider context of the classical limits of effective non-hermitian 
quantum theories. In the present letter we therefore introduce 
a generalized mean-field approximation and investigate the 
characteristic features of the dynamics resulting from the in- 
terplay of nonlinearity and non-hermiticity for a simple many- 
particle Hamiltonian of Bose-Hubbard type, describing a BEC 
in a leaking double well trap: 
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Here dj, or- are bosonic particle annihilation and creation op- 
erators for the jth mode. The onsite energies are ±8, v is 



the coupling constant and c is the strength of the onsite in- 
teraction. The additional imaginary part of the mode energy 
y describes the first mode as a resonance state with a finite 
lifetime, like, e.g., the Wannier-Stark states for a tilted optical 
lattice [9] . A direct experimental realization could be achieved 
by tunneling escape of atoms from one of the wells. Even in 
the non-hermitian case, the Hamiltonian commutes with the 
total number operator N — d\d\ + d\a2 and the number N of 
particles is conserved. The "decay" describes not a loss of 
particles but models the decay of the probability to find the 
particles in the two sites considered here. 

First theoretical results for the spectrum of the non- 
hermitian two-site Bose-Hubbard system (Q]) and a closely re- 
lated PT- symmetric system were presented in |tl(j,L!JJ]. In this 
paper we will present first results for the dynamics of this de- 
caying many-particle system with emphasis on the mean-field 
limit of large particle numbers. In order to specify the mean- 
field approximation in a controllable manner, we derive cou- 
pled equations for expectation values under the assumption 
that the system, initially in a coherent state, remains coherent 
for all times of interest. This is a direct extension of the frozen 
Gaussian approximation in flat phase space 112L l_13D to 577(2) 
coherent states, relevant to the present case as discussed be- 
low. This yields classical evolution equations for the coherent 
states parameters. 

It facilitates the analysis to rewrite the Hamiltonian (OQ) in 
terms of angular momentum operators L x = \(d\a2 -f ai^), 
L v = Y t (a\d2 - 
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iL z and cyclic permutations, as 
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The conservation of N appears as the conservation of L 
l(l + l), i.e. the rotational quantum number £ = N/2. 
The system dynamics is therefore restricted to an (N + 1)- 
dimensional subspace and can be described in terms of the 
Fock states \k,N-k), k = 0,...,Af or the 5(7(2) coherent 
states LI 4], describing a pure BEC: 
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with Xj G C. The norm, which may differ from unity, is 

(*i,*2|*i,*2) =n N , where n= |;q| 2 + |*2| 2 . 



A general discussion of the time evolution of a quantum 
system under a non-hermitian Hamiltonian 9{ — H — it with 
hermitian H and f can be found in [15]. Matrix elements 
of an operator A without explicit time-dependence satisfy the 
generalized Heisenberg equation, which in our case becomes 
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where 



is the anti-commutator. As an immediate conse- 



quence of the non-hermiticity, the norm of the quantum state 
is not conserved, h^(\\f\\\f) = -2(\|/|f |\|/) , thus the survival 
probability decays exponentially for the simple case of a con- 
stant r > 0. The time evolution of the expectation value of an 
observable (A) = (\|/|A|\|/)/(\|/|\|/) is described by the equation 
of motion 
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with the covariance A^ r = ( - [A , f ] + ) — (A) (f ) . 

For the Bose-Hubbard system these evolution equations, 
formulated in terms of the angular momentum operators with 
read (units with H = 1 are used in the following) 

±(t x ) = -2e(ly) -2c([ty,t z } + ) -2y{2A^+A^} 

i(Ly) = 2£(L x )+2c([L x ,L z } + )-2v(L z )-2j{2A 2 LyL +A 2 LytN } 



l t (L z ) = 2v(Ly) -2j{2A 2 LzLz +Al N } 
and the norm decays according to 



dt 
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In order to establish a mean-field description, we choose a 
coherent initial state \x\ ,X2) 9 i.e. a most classical state, and as- 
sume that it remains coherent for all times of interest. This as- 
sumption is, in fact, exact, if the Hamiltonian is a linear super- 
position of the generators of the dynamical symmetry group, 
i.e. for vanishing interaction c = (the proof in I114I1 can be 
directly extended to the non-hermitian case). For the inter- 
acting case c 7^ this is an approximation and the mean-field 
equations of motion are obtained by replacing the expectation 
values in the generalized Heisenberg equations of motion © 
with their values in 51/(2) coherent states ©. 

The SU(2) expectation values of the L;, / = x,y,z 9 read 
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with the abbreviations Sj = (Lj)/N for the mean values per 
particle; the expectation values of the anti-commutators fac- 
torize as 



([L / ,L / ],)=2(1-1)(L / )(L / ). 
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and ([Li,N]+) = 2N(Li). Inserting these expressions into © 
and taking the macroscopic limit N — > °° with Nc — c fixed, we 
obtain the desired non-hermitian mean-field evolution equa- 
tions: 
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These nonlinear Bloch equations are real valued and conserve 
s 2 — si + Sy + s$ = 1/4, i.e. the dynamics is regular and the 
total probability n decays as 



n = -2y(2s z + l)n. 



(11) 



Equivalently, the nonlinear Bloch equations dTUb can be 
written in terms of a non-hermitian generalization of the 
discrete nonlinear Schrodinger equation, i.e. for the time- 
evolution of the coherent state parameters x\ 9 X2- Most 
interestingly, these equations are canonical, ixj = dH/dx*-, 
ix*- = —dW/dxj, j = 1,2, where the Hamiltonian function 

is related to the expectation value of the Hamiltonian H : 
//(*!, jej,JC2>*2) = (?£) n /N and can be conveniently rewrit- 
ten in terms of the quantities \|// = e l $Xj where the (irrel- 
evant) total phase is adjusted according to (3 = — gK 2 with 

K = (|\|/i| 2 — \\\f2\ 2 )/n. The resulting discrete non-hermitian 
GPE reads 



dt V¥2 
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Similar non-hermitian mean-field equations, with the choice 



IVil 



|v|/2| , leading to different dynamics, have been 



suggested and studied before y, [5|, |_10|, LL6Q. These ad hoc 
nonlinear non-hermitian equations also appear for absorbing 
nonlinear waveguides Ill7ll . 

The dynamics of the nonlinear Bloch equations dTUb is or- 
ganized by the fixed points which are given by the real roots 





FIG. 1: (Color online) Mean-field dynamics on the Bloch sphere 
for the hermitian y = (top) and the non-hermitian case y = 0.75 
(bottom) for c = (left) and c = 2 (right) and 8 = and v = 1. 




FIG. 2: (Color online) Decay of the survival probability (full black 
curve) and the populations of site 1 (dashed red curve) and 2 (dotted 
blue curve) for an initial coherent state located at the south pole, for 
c = 0. 1, y= 0.01, v = 1 and N = 20 (left) and the relative deviations 
between many-particle and mean-field results (right). 



of the fourth order polynomial: 

4(c 2 + y 2 )4+4c£^ + (8 2 + v 2 -c 2 -y 2 )^-c8^-8 2 /4 = 0. 

(13) 
In the following we will restrict ourselves to the symmetric 
case £ = 0. Then the polynomial (fT3l) becomes biquadratic 
and the fixed points are easily found analytically. 

In parameter space we have to distinguish three different 
regions: (a) For c 2 -\-y 2 < v 2 , we have two fixed points which 
are both simple centers, (b) For |y| > |v|, we have again two 
fixed points, a sink and a source, (c) Four coexisting fixed 
points are found in the remaining region, namely a sink and 
a source (respectively two centers for y = 0), a center and a 
saddle point. Note that the index sum of these singular points 
on the Bloch sphere must be conserved under bifurcations and 
equal to two [18]. Bifurcations occur at critical parameter val- 
ues: For c 2 + y 2 = v 2 (and y ^ 0), one of the two centers (index 
+ 1) bifurcates into a saddle (index —1) and two foci (index 
+ 1), one stable (a sink) and one unstable (a source). This 
is a non-hermitian generalization of the selftrapping transi- 
tion for y = 0. The corresponding critical interaction strength 
is decreased by the non-hermiticity, i.e. the decay supports 
selftrapping . For y = ±v, the saddle (index —1) and the cen- 
ter (index +1) meet and disappear. For g = 0, we observe a 
non-generic bifurcation at y = ±v (an exceptional point 111 ill ) 
where the two centers meet and simultaneously change into a 
sink and a source. 

As an example, Fig. [T] shows the flow ([8]) on the Bloch 
sphere for v = 1 both for the hermitian y = (top) and the 
non-hermitian case y = 0.75 (bottom). For y = we observe 
the well-known selftrapping effect: In the interaction free case 
c = (upper left) we have two centers at s y = s z = 0, s x = ±^ 
and Rabi oscillations. Increasing the interaction c one of the 
centers bifurcates into a saddle (still at s z = 0) and two centers, 
which approach the poles (upper right for c = 2). The corre- 
sponding nonlinear stationary states therefore favor one of the 
wells. In the decaying system with y = 0.75 (bottom), these 
patterns are changed. For c = (lower left) we are still in 
region (a) with two centers located on the equator; however, 
they move towards s x = 0, s y = \, approaching each other. 
For c = 2 (lower right), in region (b) above the bifurcation, 



we have a center, a sink (lower hemisphere), a source (upper 
hemisphere) and a saddle. The system relaxes to a state with 
excess population in the non-decaying well, i.e. the selftrap- 
ping oscillations are damped, which is in agreement with the 
effect of decoherence in a related nonlinear two mode system 
reported in 1 19]. Finally, in region (c) only a source and a sink 
survive and the flow pattern simplifies again (not shown). The 
manifestation of the different mean-field regimes in the many 
particle system is the occurrence and unfolding of higher or- 
der exceptional points in the spectrum II ill . 

Let us finally compare the mean-field evolution with the 
full many-particle dynamics. The full quantum solution is ob- 
tained by numerically integrating the Schrodinger equation for 
the Bose-Hubbard Hamiltonian (OQ) for an initial coherent state 
with unit norm. Figure [2] shows the decay of the total survival 
probability (\|/|\|/) as a function of time for weak interaction 
(c = 0.1) and weak decay (y = 0.01) with v = 1, when ini- 
tially the non-decaying site 2 is populated. The multi-particle 
results agree with the mean-field counterpart n N on the scale 
of drawing. The deviation increases with time as can be seen 
on the right side. The probability shows a characteristic stair- 
case behavior (see also |5l [60) due to the fact that the popu- 
lation oscillates between the two sites and the decay is fast 
when site 1 is strongly populated and slow if it is empty. 
This picture is confirmed by the populations (\|/|«j^i|\|/)/Af 
and (\|/|«2<Z2|v)/N of the two sites also shown in the fig- 
ure. These quantities agree with their mean-field counterparts 
(1/2 + s z )n N /2 and (1/2 - s z )n N /2 on the scale of drawing. 
The overall decay of the norm is approximately exponential, 
37(v|/|i|/) ~ -2y/V(\|/|\|/) within region (a), as seen from (ITTb 
with sl~ = 0. 

The dynamics on the Bloch sphere in region (a) typically 
show Rabi-type oscillations. An example with parameters 
c = 0.5 and y = 0.1 is shown in Fig. [3j The mean-field os- 
cillation follows a big loop extending over the whole Bloch 
sphere. The many-particle motion oscillates with the same 
period, however, with a decreasing amplitude. This effect, 
known as breakdown of the mean-field approximation in the 
hermitian case, is due to the spreading of the quantum phase 
space density over the Bloch sphere, and can be partially cured 
by averaging over a density distribution of mean-field trajec- 
tories iEoi. 

For strong interaction, i.e. in the selftrapping region (c), we 
find an attractive fixed point, a sink, in the mean-field dynam- 
ics. An example is shown in Fig. [3] for c = 2 and y = 0.5. 
The mean-field trajectory, which started at the north pole, ap- 
proaches the fixed point at s z p = —0.433. The full many- 
particle system shows a very similar behavior. 

Further numerical investigations show that the short time 
behavior of the many-particle dynamics, as well as character- 
istic quantities such as, e.g., the half-life time, are extremely 
well captured by the mean-field description in most parameter 
ranges. 

In this letter, we have constructed a mean-field approxima- 
tion for a non-hermitian many-particle Hamiltonian, which 
can directly be generalized to other effective non-Hermitian 




FIG. 3: (Color online) Mean-field evolution of the population im- 
balance s z (t) (dashed blue curve) in comparison with the full many 
particle system for N = 20 particles (black curve) and an initial co- 
herent state located at the north pole (c = 0.5, y = 0. 1 (top) and c = 2, 
y = 0.5 (bottom) and v = 1). 



Hamiltonians. The resulting dynamics differ from the ad hoc 
non-hermitian evolution equations used in previous studies. It 
should be noted that the nonlinear Bloch equations (ITOb can 
be derived in an alternative way, based on a recently formu- 
lated number-conserving evolution equation in quantum phase 
space for M-site Bose-Hubbard systems [20], which allows 
for an immediate extension to the non-hermitian case. 
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